#include "Astronomy/Transformation.hpp"
#include <cassert>
#include <cmath>


Coordinate2D Transformation::EquatorialEcliptic(double Alpha, double Delta, double Epsilon) noexcept
{
    Coordinate2D Ecliptic;
    const double cosEpsilon = cos(Epsilon);
    const double sinEpsilon = sin(Epsilon);
    const double sinAlpha   = sin(Alpha);

    Ecliptic.X              = atan2((sinAlpha * cosEpsilon) + (tan(Delta) * sinEpsilon), cos(Alpha));
    Ecliptic.Y              = asin((sin(Delta) * cosEpsilon) - (cos(Delta) * sinEpsilon * sinAlpha));
    return Ecliptic;
}
Coordinate2D Transformation::EclipticEquatorial(double Lambda, double Beta, double Epsilon) noexcept
{
    Coordinate2D Equatorial;
    const double cosEpsilon = cos(Epsilon);
    const double sinEpsilon = sin(Epsilon);
    const double sinLambda  = sin(Lambda);

    Equatorial.X            = atan2((sinLambda * cosEpsilon) - (tan(Beta) * sinEpsilon), cos(Lambda));
    Equatorial.Y            = asin((sin(Beta) * cosEpsilon) + (cos(Beta) * sinEpsilon * sinLambda));

    return Equatorial;
}

Coordinate2D Transformation::Equatorial2Ecliptic(double Alpha, double Delta, double Epsilon) noexcept
{
    Alpha                 = HoursToRadians(Alpha);
    Delta                 = DegreesToRadians(Delta);
    Epsilon               = DegreesToRadians(Epsilon);

    Coordinate2D Ecliptic = EquatorialEcliptic(Alpha, Delta, Epsilon).RadiansToDegrees();

    if (Ecliptic.X < 0)
        Ecliptic.X += 360;

    return Ecliptic;
}

Coordinate2D Transformation::Ecliptic2Equatorial(double Lambda, double Beta, double Epsilon) noexcept
{
    Lambda                  = DegreesToRadians(Lambda);
    Beta                    = DegreesToRadians(Beta);
    Epsilon                 = DegreesToRadians(Epsilon);

    Coordinate2D Equatorial = EclipticEquatorial(Lambda, Beta, Epsilon);

    Equatorial.X            = RadiansToHours(Equatorial.X);
    Equatorial.Y            = RadiansToDegrees(Equatorial.Y);

    if (Equatorial.X < 0)
        Equatorial.X += 24;

    return Equatorial;
}
Coordinate2D Transformation::Equatorial2Horizontal(double LocalHourAngle, double Delta, double Latitude) noexcept
{
    LocalHourAngle = HoursToRadians(LocalHourAngle);

    Delta          = DegreesToRadians(Delta);
    Latitude       = DegreesToRadians(Latitude);

    Coordinate2D Horizontal;
    const double cosLatitude       = cos(Latitude);
    const double cosLocalHourAngle = cos(LocalHourAngle);
    const double sinLatitude       = sin(Latitude);

    Horizontal.X                   = RadiansToDegrees(atan2(sin(LocalHourAngle), (cosLocalHourAngle * sinLatitude) - (tan(Delta) * cosLatitude)));
    Horizontal.Y                   = RadiansToDegrees(asin((sinLatitude * sin(Delta)) + (cosLatitude * cos(Delta) * cosLocalHourAngle)));

    if (Horizontal.X < 0)
        Horizontal.X += 360;

    return Horizontal;
}

Coordinate2D Transformation::Horizontal2Equatorial(double Azimuth, double Altitude, double Latitude) noexcept
{
    // 从度转换为弧度
    Azimuth  = DegreesToRadians(Azimuth);
    Altitude = DegreesToRadians(Altitude);
    Latitude = DegreesToRadians(Latitude);

    Coordinate2D Equatorial;
    const double cosAzimuth  = cos(Azimuth);
    const double sinLatitude = sin(Latitude);
    const double cosLatitude = cos(Latitude);

    Equatorial.X             = RadiansToHours(atan2(sin(Azimuth), (cosAzimuth * sinLatitude) + (tan(Altitude) * cosLatitude)));
    Equatorial.Y             = RadiansToDegrees(asin((sinLatitude * sin(Altitude)) - (cosLatitude * cos(Altitude) * cosAzimuth)));

    if (Equatorial.X < 0)
        Equatorial.X += 24;

    return Equatorial;
}

Coordinate2D Transformation::Equatorial2Galactic(double Alpha, double Delta) noexcept
{
    Alpha = 192.25 - HoursToDegrees(Alpha);
    Alpha = DegreesToRadians(Alpha);
    Delta = DegreesToRadians(Delta);

    Coordinate2D Galactic;
    const double cosAlpha = cos(Alpha);
    const double sin274   = sin(DegreesToRadians(27.4));
    const double cos274   = cos(DegreesToRadians(27.4));

    Galactic.X            = RadiansToDegrees(atan2(sin(Alpha), (cosAlpha * sin274) - (tan(Delta) * cos274)));
    Galactic.Y            = RadiansToDegrees(asin((sin(Delta) * sin274) + (cos(Delta) * cos274 * cosAlpha)));

    Galactic.X            = 303 - Galactic.X;

    if (Galactic.X >= 360)
        Galactic.X -= 360;

    return Galactic;
}

Coordinate2D Transformation::Galactic2Equatorial(double l, double b) noexcept
{
    l -= 123;
    l = DegreesToRadians(l);
    b = DegreesToRadians(b);

    Coordinate2D Equatorial;
    const double cosl   = cos(l);
    const double sin274 = sin(DegreesToRadians(27.4));
    const double cos274 = cos(DegreesToRadians(27.4));

    Equatorial.X        = RadiansToDegrees(atan2(sin(l), (cosl * sin274) - (tan(b) * cos274)));
    Equatorial.Y        = RadiansToDegrees(asin((sin(b) * sin274) + (cos(b) * cos274 * cosl)));

    Equatorial.X += 12.25;

    if (Equatorial.X < 0)
        Equatorial.X += 360;

    Equatorial.X = DegreesToHours(Equatorial.X);

    return Equatorial;
}

double Transformation::DMSToDegrees(double Degrees, double Minutes, double Seconds, bool bPositive) noexcept
{
    // 验证参数
    if (!bPositive)
    {
        // 如果“bPositive”参数为假，则所有参数都应为非负数;
        assert(Degrees >= 0);
        assert(Minutes >= 0);
        assert(Seconds >= 0);
    }

    if (bPositive)
        return Degrees + (Minutes / 60.0) + (Seconds / 3600.0);
    else
        return -Degrees - (Minutes / 60.0) - (Seconds / 3600.0);
}

Coordinate3D Transformation::SphericalToCartesian(const Coordinate3D &spherical) noexcept
{
    const double Phi      = spherical.X;
    const double Theta    = spherical.Y;
    const double R        = spherical.Z;

    const double cosTheta = cos(Theta);

    Coordinate3D cartesian{0};
    cartesian.X = R * cosTheta * cos(Phi);
    cartesian.Y = R * cosTheta * sin(Phi);
    cartesian.Z = R * sin(Theta);

    return cartesian;
}

Coordinate3D Transformation::CartesianToSpherical(const Coordinate3D &cartesian) noexcept
{
    const double x = cartesian.X;
    const double y = cartesian.Y;
    const double z = cartesian.Z;

    // Phi = X,Theta = Y,R = Z;
    Coordinate3D spherical{0};

    spherical.Z = sqrt(x * x + y * y + z * z);
    spherical.Y = asin(z / spherical.Z);
    spherical.X = atan2(y, x);
    return spherical;
}